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ABSTRACT 


There  are  several  perennial  issues  concerning  reliability  estimation  in  Operational 
Testing  and  Evaluation  (OT&E).  Two  of  these  include,  "how  should  one  model  the 
underlying  failure  distribution  of  a  continuous-time  system,"  and,  "how  can  a  testing 
agency  use  information  from  DT  in  order  to  reduce  OT  resource  requirements.” 

In  the  former  issue,  some  OT&E  analysts  have  questioned  whether  or  not  the 
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First,  this  thesis  compares  the  exponential  and  Weibull  failure  distributions  in  terms  of 
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MTTF  of  a  system  is  at  least  as  large  as  the  minimum  acceptable  value,  and  also  in  terms 
of  the  actual  confidence  level  associated  with  the  lower  confidence  level  procedure  when 
the  system  has  an  increasing  (or  decreasing)  failure  rate  function.  Second,  the  thesis 
examines  the  behavior  of  an  estimator  for  the  relationship  between  DT  and  OT  failure  data 
using  a  Monte  Carlo  simulation.  Finally,  the  thesis  introduces  a  hierarchical  Bayes 
approach  for  the  estimation  of  the  relationship  between  DT  and  OT  failure  data  when  a 
gamma  prior  distribution  is  assumed. 
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EXECUTIVE  SUMMARY 


There  are  several  perennial  issues  concerning  reliability  estimation  in  Operational 
Testing  and  Evaluation  (OT&E).  One  of  these  is,  "how  should  one  model  the  underlying 
failure  distribution  of  a  continuous-time  system,"  and  another  is,  "how  can  a  testing 
agency  use  information  from  Developmental  Testing  (DT)  in  order  to  reduce  Operational 
Testing  (OT)  resource  requirements." 

The  purpose  of  this  thesis  is  to  (1)  investigate  the  behavior  of  a  statistical  procedure 
used  to  determine  the  length  of  time  needed  to  test  a  continuously  operating  system  when 
the  system  under  test  is  assumed  to  have  an  exponential  time  to  failure  distribution,  but 
when  in  fact  the  system  does  not  have  an  exponential  time  to  failure  distribution  and  rather 
may  have  a  time  to  failure  distribution  that  more  closely  resembles  a  Weibull  distribution; 
(2)  to  examine  the  behavior  of  a  maximum  likelihood  estimator  for  the  relationship 
between  DT  and  OT  failure  rates  using  a  Monte  Carlo  simulation;  and  (3)  to  introduce  a 
hierarchical  Bayes  method  for  estimation  of  the  relationship  between  the  DT  and  OT 
failure  rates.  In  the  first  part  of  the  thesis,  the  parameters  of  the  procedure  are  chosen 
based  on  the  assumption  of  an  exponential  time  to  failure  distribution.  The  behavior  of  the 
procedure  is  then  investigated  under  the  assumption  that  the  time  between  failures  are 
independent,  having  a  Weibull  distribution  with  the  same  mean  time  to  failure.  The 
intention  is  to  test  whether  or  not  use  of  the  Weibull  distribution  will  result  in  fewer  test 
hours  to  verify  that  a  system  under  test  has  a  mean  time  to  failure  which  meets  or  exceeds 
a  minimum  acceptable  value  at  a  given  level  of  confidence,  and  also  to  test  whether  or  not 

I 

the  OTAs  are  misstating  the  confidence  they  report  in  their  operational  testing  results 
when  the  exponential  distribution  is  used  to  model  a  system  whose  true  failure  rate 
function  is  that  of  a  Weibull  distribution.  In  the  second  part  of  the  thesis,  a  maximum 


likelihood  estimator  is  examined  in  terms  of  the  consistency  and  variability  of  the  OT 
failure  rate  derived  using  that  estimator  compared  to  other  methods  of  establishing  the  OT 
failure  rate.  Finally,  in  the  third  part  of  the  thesis,  a  hierarchical  Bayes  model  is  proposed 
for  the  ultimate  purpose  of  being  able  to  describe  the  mean  and  variance  associated  with  a 
random  relationship  between  DT  and  OT  failure  rates. 

The  comparison  of  the  exponential  and  Weibull  failure  models  suggests  that 
operational  testing  agencies  can,  under  certain  circumstances,  reduce  the  amount  of  time 
needed  to  demonstrate  that  the  mean  time  to  failure  (MTTF)  of  a  continuously  operating 
system  is  at  least  as  large  as  the  minimum  acceptable  value  (MAV)  if  the  system's  true 
underlying  failure  distribution  is  Weibull  with  shape  parameter  |3  *  1.  However,  the  shape 
parameter  is  seldom,  if  ever,  known  or  even  estimable  without  expending  valuable  test 
resources,  so  it  is  unlikely  that  an  operational  testing  agency  (OTA)  could  actually  save 
any  time  or  resources  with  this  methodology.  Moreover,  from  the  perspective  of 
confidence  levels,  the  analysis  suggests  that  confidence  in  one's  result  is  not  sacrificed 
when  the  exponential  distribution  is  used  to  model  a  system  which  has  a  non-constant 
failure  rate,  since  it  appears  that  use  of  the  exponential  distribution  on  a  truly  Weibull 
system  seems  to  cause  the  OTA  to  overestimate  the  number  of  hours  needed  to  test  the 
MAV  for  the  MTTF  to  a  given  confidence  level,  and  also  to  underestimate  the  confidence 
it  reports  in  its  results  under  the  incorrect  assumption  of  an  exponential  distribution.  The 
thesis  also  points  out,  however,  that  these  results  may  not  be  true  in  general.  Therefore,  it 
is  concluded  that  there  is  currently  little  to  be  gained  by  attempting  to  model  the  time 
between  system  failures  with  the  Weibull  distribution  as  previously  suggested,  particularly 
if  very  few  observations  of  time  to  failure  arc  available. 

The  second  analysis  indicates  that  the  maximum  likelihood  method  for  estimating  a 
system's  OT  failure  rate  given  that  system's  DT  failure  rate  as  well  as  failure  data  on 


previous  similar  systems  produces  a  reasonable  estimate  of  the  actual  OT  failure  rate 
(when  the  OT  failure  rate  behaves  as  expected).  In  fact,  the  estimate  provided  by  the 
maximum  likelihood  method  appears  to  be  nearly  as  good  as  direct  observation  of  the  OT 
failure  rate  obtained  through  actual  testing.  Thus,  if  it  can  be  shown  that  inferences  about 
the  failure  rate  of  a  system  can  in  fact  be  made  on  the  basis  of  previous  similar  systems, 
then  the  maximum  likelihood  method  of  predicting  a  new  system's  OT  failure  rate  from 
previous  DT  failure  data  has  some  potential  utility.  In  particular,  this  OT  failure  rate 
prediction  could  be  used  in  ways  which  would  aid  decisions  as  to  when  to  begin  OT  as 
well  as  to  augment  follow-on  OT  in  an  effort  to  reduce  needed  testing  resources. 

The  third  analysis  is  incomplete.  It  was  hoped  that  the  hierarchical  Bayes  method 
given  would  provide  a  means  by  which  to  estimate  the  relationship  between  the  DT  and 
OT  failure  rates  when  the  relationship  is  assumed  to  be  a  random  variable  with  some  prior 
distribution.  Regretfully,  the  equations  derived  for  the  maximum  likelihood  estimators  for 
the  parameters  of  the  gamma  prior  have  proved  to  be  difficult  to  solve  satisfactorily.  If  a 
solution  to  these  equations  were  available,  it  could  be  used  to  predict  future  failure  rates 
by  means  of  a  Bayesian  calculation.  Investigation  of  such  a  procedure  is  left  for  a  later 
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I.  INTRODUCTION 


There  are  several  perennial  issues  concerning  reliability  estimation  in  Operational  Test 
and  Evaluation  (OT&E).  In  continuous  systems,  such  as  radios,  radars  and  other  items 
that  operate  continuously  while  in  use,  one  perennial  issue  is  "how  should  one  model  the 
underlying  failure  distribution  of  the  system?"  Another  is  "what  can  be  done  to  reduce  the 
amount  of  testing  resources  needed  during  OT&E?" 

Historically,  OT&E  organizations,  such  as  Commander,  Operational  Test  and 
Evaluation  Force  (COTF),  have  used  the  exponential  distribution  to  model  the  behavior  of 
the  reliability  of  a  continuous-time  type  system.  The  assumption  of  an  exponential  time  to 
failure  distribution  results  in  a  required  test  period  of  three  times  a  given  minimum 
acceptable  value  (MAV)  for  the  80  percent  lower  confidence  limit  for  the  mean  time  to 
failure  (MTTF).  As  long  as  fewer  than  two  failures  occur  during  this  time,  then  COTF 
will  report  80  percent  confidence  that  the  true  MTTF  of  the  system  is  at  least  as  large  as 
the  MAV  (COTF,  para.  204.a.(3),  1992).  COTF  analysts  have  questioned  whether  or  not 
the  exponential  distribution  is  the  best  underlying  failure  model  to  use,  and  have 
considered  investigating  alternative  distributions  (Madson,  1993).  One  possible 
alternative  in  question  is  the  Weibull  distribution,  which  contains  the  exponential 
distribution  as  a  special  case.  The  main  difference  between  the  two  underlying  failure 
distributions  is  that  the  exponential  distribution  has  a  constant  failure  rate  function 
throughout  the  lifetime  of  the  system,  and  the  Weibull  distribution  has  a  changing  failure 
rate  function  (either  increasing,  constant  or  decreasing)  during  the  lifetime  of  the  system. 
The  idea  is  that  the  Weibull  distribution  might  provide  a  "better  fitting"  model  for  certain 
systems  which  might  then  reduce  the  amount  of  test  time  needed  to  demonstrate,  to  a 
given  confidence  level,  that  the  true  MTTF  of  a  system  is  at  least  as  large  as  the  MAV. 


A  third  perennial  issue  is  that  of  combining  Developmental  Testing  (DT)  failure  data 
and  Operational  Testing  (OT)  failure  data  in  a  way  that  would  reduce  the  overall  amount 
of  resources  needed  to  establish  the  reliability  of  a  system.  For  example,  Gaver  and 
Jacobs  (1994,  p.  1)  describe  the  following  problem  setting: 

A  test  and  evaluation  analyst  has  in  his  possession  historical  (time-between-failure 
or  equivalent)  data  for  both  the.. .DT... and.. .OT... phases  of  a  number  of  projects  that 
are  roughly  comparable:  they  are  all  C3/I  systems,  perhaps.  It  is  believed  that  there 
may  be  a  usefully  exploitable  relationship  between  the  DT  and  OT  data.  If  this 
relationship  could  be  quantified  then  perhaps  it  could  be  used  to  augment  DT  data 
for  a  new  project  of  the  same  type,  thereby  obtaining  some  anticipation  of  the  OT 
data  for  that  new  project... 

Gaver  and  Jacobs  (1994,  p.  5)  go  on  to  devise  an  estimator  for  the  relationship 
between  the  DT  data  and  OT  data  and  describe  its  application  thusly: 

The  ultimate  application  of  the  above  (relationship)  is  to  (a)  anticipate  failure 
patterns  during  Operational  Testing,  given  data  on  failures  during  Developmental 
Testing;  this  might  aid  in  the  decision  as  to  when  to  begin  OT;  and  (b)  to 
strengthen,  and  reduce  uncertainty  of,  the  post-DT  estimates  of...(system 
reliability)...by  incorporating  the  DT  data. 

Of  course,  in  the  preceding  example,  it  is  assumed  that  the  relationship  between  the 
DT  and  OT  failure  rates  is  constant  across  all  "similar"  systems.  It  is  far  more  likely  that 
the  relationship  follows  some  probability  distribution.  Thus,  a  hierarchical  approach  for 
the  description  of  this  relationship  may  be  appropriate.  One  such  approach  is  outlined  in 
the  last  part  of  this  thesis. 

The  notion  of  combining  DT  and  OT  data  has  historically  been  an  anathema  to  those 
involved  in  OT&E.  The  primary  reason  for  this  is  that  the  Operational  Testing  Agencies 
(OTAs)  have,  rightly  so,  desired  to  be  completely  independent  of  the  Developmental 
Agencies  (DAs).  However,  with  the  tight  budgets  facing  all  military  organizations  of  late, 
it  may  be  time  to  re-visit  the  idea  of  using  previous  DT  data  to  supplement  OT  data  in  an 
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effort  to  reduce  the  amount  of  testing  resources  needed  to  establish  a  system's  reliability  in 
OT&E 

The  first  goal  of  this  thesis  is  to  compare,  through  a  sensitivity  analysis,  the 
exponential  and  Weibull  failure  distributions  in  terms  of  the  amount  of  test  time  needed  to 
demonstrate,  to  a  given  level  of  confidence,  that  the  true  MTTF  of  a  system  is  at  least  as 
large  as  the  MAV,  and  to  find  the  actual  confidence  level  associated  with  the  lower 
confidence  limit  proced’  e  when  the  system  has  an  increasing  (or  decreasing)  failure  rate 
function.  The  former  tests  whether  or  not  use  of  the  Weibull  distribution  will  result  in 
fewer  test  hours  to  establish  reliability  with  a  given  confidence  level.  The  latter  tests 
whether  or  not  the  OTA  is  overstating  (or  understating)  the  confidence  it  reports  in  its 
results  when  the  exponential  distribution  is  assumed  but  the  system's  true  failure  rate 
follows  a  Weibull  distribution. 

The  second  goal  of  this  thesis  is  to  examine  the  behavior  of  the  estimator  for  the 
relationship  between  DT  and  OT  data  devised  by  Gaver  and  Jacobs  (1994)  using  a  Monte 
Carlo  simulation.  The  behavior  to  be  examined  includes  the  consistency  and  variability  of 
the  OT  failure  rate  derived  using  this  estimator  compared  to  those  of  other  methods  of 
establishing  the  OT  failure  rate,  including  direct  observation  of  the  OT  data  and  a 
weighted  average  between  direct  observation  and  the  Gaver-Jacobs  estimate  of  the  OT 
failure  rate. 

The  third  goal  of  this  thesis  is  to  introduce  a  hierarchical  Bayesian  approach  for  the 
estimation  of  the  relationship  between  DT  and  OT  data  using  a  gamma  prior  distribution. 

Note  that  the  methods  used  for  borrowing  DT  information  to  assist  in  predicting  OT 
information,  or  actual  later  field  experience,  can  also  be  used  in  other  situations.  For 
example,  suppose  an  upgrade  of  an  existing  system  is  to  be  evaluated.  Then  actual  data 
from  the  existing  system  may  be  scaled  to  the  new  system  in  a  manner  analogous  to  that 


already  described  for  DT-to-OT  scaling.  Note  also  that  the  general  procedure  described 
can  be  adapted  to  other  operational  parameters.  An  example  might  be  an  effectiveness 
parameter  such  as  projectile  dispersion  around  an  aim  point,  technically  speaking,  a 
variance  parameter.  Under  current  conditions  of  declining  resources  such  an  approach 
could  be  of  considerable  value  in  enhancing  the  accuracy  of  many  types  of  test  data. 
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n.  COMPARISON  OF  EXPONENTIAL  AND  WEIBULL  DISTRIBUTIONS 

A.  EXPONENTIAL  MODEL 

One  of  the  simplest  and  most  popular  ways  to  model  the  distribution  of  the  random 
time  to  failure  of  a  continuous  time  system  is  with  the  exponential  distribution.  Among 
the  reasons  for  the  popularity  of  the  exponential  distribution  are  its  simplicity,  ease  of  use, 
and  the  memoryless  property  which  says  that  the  probability  of  failure  within  the  next  t 
time  units  is  independent  of  the  frequency  of  any  previous  failures. 

Let  T  =  random  time  to  failure  of  a  system, 
t  =  length  of  time  system  is  in  operation, 

7  =  failure  rate. 

The  time  to  failure,  T,  of  a  system  follows  an  exponential  failure  distribution  if  the 
probability  density  function  (PDF)  is 


f(t)  =  -ye-^t  >0 

and  if  its  cumulative  distribution  function  (CDF)  is 


(2.1) 


F(t)  =  l-e~7‘,t>0 
=  0,t<0 


(2.2) 


If  the  times  between  system  failures  are  independent  and  have  the  exponential 
distribution  specified  above,  then  the  system's  mean  time  to  failure  (MTTF)  is 

MTTF  =  E[T]  =  -  (2.3) 

T 

and  its  failure  rate  function  is 


h(t)  = 


f(t) 

l-F(t) 


=  7. 


(2.4) 
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This  distribution  is  used  by  the  analysts  and  Operational  Test  Directors  (OTDs)  at 
COTF  primarily  to  determine  the  length  of  time  needed  to  test  a  given  continuous  time 
system  in  order  to  demonstrate  that  the  true  MTTF  of  the  system  is  greater  than  or  equal 
to  the  MAV  for  the  MTTF  with  approximately  known  confidence.  The  basic  assumption 
associated  with  the  exponential  distribution  is  that  the  system  under  scrutiny  exhibits 
steady  state  performance;  that  is,  the  time  between  failures  are  independent  and  have  a 
common  exponential  distribution,  in  other  words,  it  is  assumed  that  the  failure  rate 
function,  h(t),  is  constant  over  time:  there  is  no  "reliability  growth,"  or  degradation. 

COTF's  usual  policy  is  to  conduct  sequential  tests.  That  is,  they  will  test  one  system 
for  a  pre-determined  length  of  time.  If  it  fails  before  the  time  on  test  expires  then  they  will 
repair  the  system  and  continue  the  test.  If  the  repaired  system  fails  before  the  remaining 
time  on  test  then  the  system  fails  the  test,  else  it  passes.  In  other  words,  if  there  are  zero 
or  one  failures  during  the  time  on  test,  then  the  system  passes  (COTF,  para.  602.a,  1992). 

Let  a  =  confidence  level, 

0  =  MAV  for  the  MTTF, 

T2  =  time  to  2nd  failure, 

N(t)  =  number  of  failures  to  occur  in  [0,  tj. 

The  time,  t,  needed  for  100(a)  percent  confidence,  with  no  more  than  one  failure 
allowed  during  testing,  required  to  demonstrate  that  the  system  meets  the  MAV  for  the 
MTTF,  0,  is  determined  as  follows: 


1  -  a  >  P(T2  >  t)  =  P(N(t)  <  1) 
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Standard  procedure  at  COTF  is  to  set  the  confidence  level  at  HO  percent  and  to  allow 
for  a  maximum  of  one  failure  to  occur  during  the  operaaonai  test.  Soluuon  of  equauon 
(2.5)  with  a  =  0.80  gives  a  time  to  test  t  =  39. 1  In  other  w^rds.  COTF  will  plan  to  test  a 
continuous  time  system  for  a  period  of  time  which  is  three  times  greater  than  the  MAV  for 
the  MTTF.  and  if  the  system  has  fewer  than  two  failures  in  this  time,  then  the  system 
passes  the  test  at  the  80  percent  confidence  level.  Note  that  other  combinations  of  test 
times  and  acceptable  numbers  of  failures  can  yield  the  same  confidence,  but  the  above  has 
been  chosen  as  a  practical  compromise. 

For  example,  if  the  MAV  for  the  MTTF  is  1000  hours,  COTF  will  test  the  system  for 
3000  hours  and  allow  at  most  one  failure  to  occur.  If  the  system  passes  this  test,  COTF 
will  report  that  it  is  80  percent  confident  that  the  system's  true  MTTF  meets  or  exceeds 
the  MAV. 

Implicit  in  the  preceding  model,  however,  is  the  notion  that  if  the  expected  value  of 
the  time  between  failures  for  one  system  is  greater  than  or  equal  to  the  expected  value  of 
the  time  between  failures  for  a  second  system,  then  the  probability  of  having  J  or  fewer 
failures  in  a  renewal  process  observed  for  a  fixed  time  for  the  first  system  is  at  least  as 
large  as  that  for  the  second  system  given  the  same  fixed  time. 

Let  X  =  time  to  1st  failure  of  system  1,  X  -  exp(7,), 

Y  =  time  to  1st  failure  of  system  2,  Y  -  expCy2),  yj  £  72- 

Proposition  1:  If  E[X]  >  E[Y]  then  P(Nx(t)  Sl)>  P(NY(t)  <  1)  V  t. 

Proof:  y,  <  y2  =»  E[X]  =  1/y,  >  E[Y]  =  l/y2 

7,  <  72  =t>  Fx(t)  =  e-7,t  >  Fy(t)  =  e_T2t  Vt>0 


1  An  interesting,  and  simpler,  alternative  method  for  determining  t  exists,  involving  the  chi-square 
distribution.  See  Appendix  A  for  details  of  the  alternative  method. 
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Let  Nx(t)  =  the  number  of  renewals  that  occur  m  a  renewal 
process  having  inter-arrival  distribution  Fx(t), 

N y ( t )  =  the  number  of  renewals  that  occur  in  a  renewal 
process  having  inter-arrival  distribution  FY(t). 

then  Fx(t)>  Fy(t)  =»  Nx(t)  ^Nyft)  V  t  >0  (Ross,  p.257,  1983), 

which  means  that  P(Nx(t)  <  1)  >  P(NY(t)  <  1).  Q.ED. 

Note  that  Proposition  1  may  not  be  true  for  all  eligible  distributions.  However,  tht 
fact  that  it  is  true  for  the  exponential  distribution  enables  one  to  use  the  COTF  procedure. 

B.  ON  CONFIDENCE  LEVELS 

One  might  ask  at  this  point  how  COTF  defines  "level  of  confidence."  As  used  in 
equation  (2.5),  the  confidence  level  is  actually  the  "producer's  risk."  When  COTF  desires 
to  test  a  system's  MTTF  with  80  percent  confidence,  this  means  they  are  willing  to  assume 
a  20  percent  (or  smaller)  risk  of  incorrectly  recommending  a  system  for  production  based 
on  the  results  of  their  test 

Let  Nx(t)  and  NY(t)  denote,  respectively,  the  number  of  renewals  to  occur  in  two 
renewal  processes  having  distribution  Fx(t)  and  FY(t).  Ross  (p.  257,  1983)  shows  that  if 
Fx(t)  >  Fy(t)  then  P(Nx(t)  <  1)  >  P(NY(t)  £  1).  Since  Fx(t)  >  Fv(t)  implies  that  E[X]  > 
E[Y],  the  result  shows  that  under  this  additional  assumption,  if  the  true  MTTF  of  the 
system  is  less  than  the  MAV,  8,  and  the  desired  level  of  confidence  is  a,  then  Ross's  result 
implies  that  P(N(t)  <  1)  <  1  -  a  and  P(N(t)  >  1)  >  a. 

In  other  words,  if  COTF  conducted  the  same  test  many  times  on  a  system  with 
exponential  times  between  failures,  and  if  the  true  MTTF  were  actually  less  than  the 
MAV,  8,  then  the  observed  time  for  the  second  failure  to  occur,  T2,  would  fall  below  38 
more  than  80  percent  of  the  time  (in  which  case  COTF  would  correctly  reject  the  system). 
Consequently,  less  than  20  percent  of  the  time  the  observed  time  for  the  second  failure  to 
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occur,  T2,  would  lie  above  39  (in  which  case  they  would  incorrectly  accept  the  system). 
Symbolically,  COTF  desires: 


P(accept  the  system  I  true  MTTF  <  MAV)  <,  .20. 

There  is  an  interesting  implication  to  this  requirement.  In  order  to  have  any 
appreciable  chance  of  passing  a  test  set  up  according  to  COTF  standards,  the  DA  must 
build  its  system  such  that  its  true  MTTF  is  much  greater  than  the  MAV  (Keller,  1993).  If 
the  system's  MTTF  is  equal  to  the  MAV,  then  it  has  only  a  20  percent  chance  of  passing 
the  test. 

Note  carefully:  the  above  does  uqi  say  that  if  the  test  leads  to  acceptance,  i.e.,  that 
fewer  than  two  failures  are  observed  in  three  times  the  MAV,  then  the  probability  is  80 
percent  that  the  true  MTTF  is  greater  than  the  MAV.  This  statement  is  specifically  not 
allowed  in  the  present  framework. 

C.  WEIBULL  MODEL 

! 

An  alternative,  and  more  general,  way  to  model  the  probability  of  failure  of  a 
continuous  system  is  with  the  Weibull  distribution.  The  Weibull  distribution  is  useful 
when  the  system  under  scrutiny  displays  either  an  increasing  or  decreasing  failure  rate  over 
i  time.  An  example  of  a  system  with  an  increasing  failure  rate  function  would  be  one  that 

"wears  out"  over  time.  A  system  that  exhibits  "wear  out"  would  tend  to  have  a  higher 
failure  rate  in  old  age  than  an  identical  system  which  is  younger.  An  example  of  a  system 
with  a  decreasing  failure  rate  function  would  be  one  that  "wears  in"  over  time.  Certain 
!  electronic  devices  may  exhibit  "wear  in"  tendencies,  at  least  initially. 

I 
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Let  T  -  random  time  to  failure  of  a  system, 
t  =  length  of  time  system  is  in  operation, 

X  =  scale  parameter  (characteristic  life), 

(3  =  shape  parameter. 

The  time  to  failure,  T,  of  a  system  follows  a  Weibull  failure  distribution  if  the 
probability  density  function  (PDF)  is 

f  (t)  =  pxV-'e-^’ ,  t  >  0,  p  >  OA  >  0  (2.6) 

and  if  its  cumulative  distribution  function  (CDF)  is 


F(t)  =  l-e^,t£0,P>0,X>0 
=  0 .  t  <  0. 


(2.7) 


If  the  times  between  system  failures  are  independent  and  have  the  Weibull  distribution 
specified  above,  then  the  mean  time  to  failure  (MTTF)  of  the  system  is 


MTTF  =  E[T] 


and  its  failure  rate  function  is 


h(t)  = 


f(t) 

l-F(t) 


P*ptM. 


(2.8) 


(2.9) 


This  distribution  has  been  considered  as  a  replacement  for  the  exponential  distribution 
with  the  hope  that  it  might  more  precisely  model  many  of  the  systems  under  test  by  COTF. 
In  this  model,  X  is  analogous  to  the  exponential  failure  rate,  y.  The  shape  parameter,  (3, 
determines  whether  the  system's  failure  rate  is  increasing,  decreasing,  or  constant  over 
time.  When  {3  is  less  than  one,  the  failure  rate  function  is  decreasing,  when  (3  is  greater 
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than  one,  the  failure  rate  function  is  increasing,  and  when  p  is  equal  to  one,  the  failure  rate 
function  is  constant.  Note  that  when  P  is  equal  to  one,  the  Weibull  distribution  is  the 
exponential  distribution. 

There  is  a  useful  relationship  between  the  exponential  M  i  ll-,  8  =  1/y,  and  the 
Weibull  parameters  A  and  P  (when  the  exponential  and  Weibull  MTTFs  are  equal): 


This  relationship  will  be  useful  when  comparing  the  performance  of  the  two  distributions. 


Let  a  =  confidence  level, 

T2  =  time  to  2nd  failure, 

N(t)  =  number  of  failures  to  occur  in  [0,  t]. 

The  time,  t,  needed  for  100(a)  percent  confidence,  with  no  more  than  one  failure 
allowed  during  testing,  required  to  demonstrate  that  the  system  meets  the  MAV  for  the 
MTTF,  6,  may  be  determined,  using  an  adaptation  of  the  test  for  the  exponential 
distribution,  equation  (2.5),  to  the  Weibull  distribution  as  follows: 


1  -  a  >  P(T2  >  t)  =  P(N(t)  £1)  =  P(N(t)  =  0)  +  P(N(t)  =  I) 
P(N(t)  =  0)  =  e_<Xl>> 

l  t 

P(N(t)  =  1)  =  J  FTj  (t  -  x)fTj  (x)dx  =  pAp  J 


l-a>e-(Xj)°  + 


(2.i  n 


Solution  of  equation  (2.11),  with  a,  (3,  and  X  specified  gives  the  appropriate  time  to 
test,  t.  Note,  however,  that  there  is  no  closed  form  solution  for  t  in  this  inequality.  It 
must  be  solved  by  numerically  integrating  for  a  number  of  t- values. 

Implicit  in  the  preceding  model,  again,  is  the  notion  that  if  the  expected  value  of  the 
time  between  failures  for  one  system  is  greater  than  or  equal  to  the  expected  value  of  the 
time  between  failures  for  a  second  system,  then  the  probability  of  having  J  or  fewer 
failures  in  a  renewal  process  observed  for  a  fixed  time  for  the  first  system  is  at  least  as 
large  as  that  for  the  second  system  given  the  same  fixed  time. 

Let  X  =  time  to  1st  failure  of  system  1,  X  -  Weibull(  P,  X,), 

Y  =  time  to  1st  failure  of  system  2,  Y  -  Weibull(  (3,  Xj),  X,  £  Xj ,  p  >  0. 

Proposition  2:  If  E[X]  >  E[Y]  then  P(Nx(t)  £  1)  >  P(NY(t)  <  1)  V  t. 

Proof:  Xj  S  Xj  => 

Xj  £  X2  => 

Let  Nx(t)  =  the  number  of  renewals  that  occur  in  a  renewal 
process  having  inter  arrival  distribution  Fx(t), 

NY(t)  =  the  number  of  renewals  that  occur  in  a  renewal 
process  having  inter  arrival  distribution  FY(t). 

then  Fx  (t)  >  Fy  (t)  =>  Nx(t)  ^  NY(t)  V  t  >  0  (Ross,  p.257, 1983) 

which  means  that  P(Nx(t)  <  1)  >  P(NY(t)  <  1).  Q.ED. 

D.  TEST  TIME  NEEDED  TO  DEMONSTRATE  A  THRESHOLD 

Suppose  COTF  is  testing  a  new  radar  system  whose  MAV  for  MTTF  is  1000  hours. 
Suppose  further,  for  the  time  being,  that  the  radar  system's  true  MTTF  is  in  fact  1000 
hours,  i.e„  y  =  0.001.  Using  the  assumption  of  a  constant  failure  rate  function,  and  the 
thumb  rule  of  80  percent  confidence  with  no  more  than  one  failure,  COTF  would  then  test 


E[X,-TT^  m]'lor 
px,  px2 

Fxd^e-^  2:  FrO^e-^  Vt20,p20 
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the  system  for  3000  hours  and  allow  at  most  one  failure  to  occur.  If  the  system  passes  this 
test,  COTF  will  report  that,  with  80  percent  confidence,  the  radar  system  MTTF  meets  or 
exceeds  the  MAV. 

Now  suppose  COTF  wishes  to  model  this  system's  behavior  with  a  Weibull 
distribution.  Several  choices  for  the  shape  parameter,  P,  are  available.  Assume  that, 
depending  on  the  platform  on  which  the  radar  system  will  be  deployed,  the  shape 
parameter  may  vary  anywhere  between  0.5  and  3.0.  The  relationship  derived  in  equation 
(2.10)  can  now  be  used  to  match  the  MTTF  of  the  exponential  distribution,  1/7,  to  the 
MTTF  of  a  Weibull  distribution  given  various  shape  parameters,  P,  as  seen  in  Table  1. 


TABLE  1 

RELATIONSHIP  OF  EXPONENTIAL  FAILURE  RATE,  7,  TO  WEIBULL  SCALE 
PARAMETER,  X,  TO  MATCH  MTTFs  FOR  DIFFERENT  WEIBULL  SHAPE 


PARAMETERS,  p 


Exponential 
failure  rate,  7 

Weibull 

shape  parameter,  P 

m/p) 

Weibull 

scale  parameter,  A. 

0.001 

0.5 

1 

0.002 

0.001 

1 

1 

0.001  | 

0.001 

2 

1.7725 

0.000886 

0.001 

3 

2.6789 

0.000893 

Now  it  is  possible  to  use  equation  (2.1 1)  to  derive  the  appropriate  test  times  for  the 
preceding  values  of  p  as  shown  in  Table  2,  again  assuming  a  confidence  level  of  80 
percent,  and  allowing  for  at  most  one  failure  to  occur.2 


2The  author  solved  equation  (2.1 1)  using  MathCad  version  3.1,  varying  t  in  one  hour  increments  until  the 
inequality  was  satisfied. 
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TABLE  2 

TEST  TIMES  ASSOCIATED  WITH  VARIOUS  WEEBULL 
_  SHAPE  PARAMETERS,  3 _ 


p 

X 

t  (hrs) 

0.5 

0.002 

2950 

1 

0.001 

3000 

2 

0.000886 

2613 

3 

0.000893 

2435 

This  analysis  suggests  that,  everything  else  being  equal,  the  assumption  of  an 
exponential  failure  distribution  (P  =  1)  generates  the  longest  test  time  to  demonstrate  that 
the  system  meets  or  exceeds  the  MAV  for  MTTF  at  an  80  percent  confidence  level.  In 
other  words,  if  it  is  reasonable  to  assume  that  the  system  under  test  follows  a  Weibull  time 
to  failure  distribution,  and  if  the  OTA  can  identify  a  reasonable  {$  value  for  the  system 
(where  p  *  1),  then  by  modeling  the  system  with  a  Weibull  failure  distribution  the  OTA 
may  be  able  to  reduce  the  required  number  of  hours  needed  to  demonstrate  that  the 
system  meets  or  exceeds  the  MAV  for  MTTF  for  a  Weibull  distribution  with  the  same 
shape  parameter  while  maintaining  the  same  level  of  confidence  in  its  result 

E.  USE  OF  EXPONENTIAL  MODEL  ON  A  WEIBULL  PROBLEM 

Now  suppose  COTF  uses  its  traditional  rule  for  determining  time  needed  to  test  the 
radar  system.  Earlier  it  was  shown  that  for  an  exponential  failure  distribution  with  a  MAV 
for  MTTF  of  1000  hours,  it  was  necessary  to  test  the  radar  system  for  3000  hours, 
allowing  for  at  most  one  failure,  to  obtain  80  percent  confidence  that  the  true  MTTF  met 
or  exceeded  the  MAV.  Suppose  it  is  later  discovered  that  the  true  underlying  failure 
distribution  was  actually  Weibull  with  shape  parameter  P  *  1.  What  can  then  be  said 
about  the  actual  level  of  confidence  achieved  from  the  test?  Is  it  still  80  percent? 

One  way  to  explore  this  is  to  hold  the  test  time  constant  at  3000  hours,  and  determine 
the  probability  of  passing  the  test  (  1  -  a  )  and  also  the  level  of  confidence  (  a  )  achieved 
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when  the  underlying  failure  distribution  is  Weibull.  This  can  be  accomplished  using 
equation  (2.1 1)  by  fixing  the  time  to  test,  t,  at  3000  hours,  specifying  P  and  X,  and  then 
solving  for  a.  Results  are  shown  in  Table  3. 3 

This  analysis  seems  to  suggest  that  no  confidence  is  lost  by  assuming  an  exponential 
failure  distribution  when  the  true  underlying  failure  distribution  is  Weibull  with  parameter 
P  >  1.  This  makes  intuitive  sense.  As  was  suggested  back  in  Table  2,  for  any  P  >  1,  the 
required  test  time  is  reduced,  thus  if  the  system  is  tested  longer  than  necessary  to  achieve 
80  percent  confidence,  as  occurs  in  Table  3  when  P  >  1,  then  it  naturally  follows  that  the 
actual  confidence  level  in  the  result  of  the  test  should  increase. 


TABLE  3 

ACTUAL  CONFIDENCE  LEVELS  ACHIEVED  AS  A  FUNCTION  OF 
THE  SHAPE  PARAMETER,  p,  WHEN  EXPONENTIAL 


DISTRIBUTION  IS  USED  TO  MODEL  A  WEIBULL  (3,X)  PI 

ROBLEM 

P 

X 

t  (hrs) 

P(pass  test) 

(1-cc) 

Conf.  level 
(a) 

0.5 

0.002 

3000 

0.196 

0.804 

1 

0.001 

3000 

0.200 

0.800 

2 

0.000886 

3000 

0.097 

0.903 

3 

0.000893 

3000 

0.028 

0.972 

Unfortunately,  the  results  of  this  analysis  may  not  hold  under  all  circumstances.  The 
reason  for  this  is  because  Propositions  1  and  2  no  longer  apply.  In  other  words,  the  fact 
that  E[X]  >  E[Y]  no  longer  implies  that  P(Nx(t)  <  1)  >  P(NY(t)  <  1). 

Let  X  =  time  to  1st  failure  of  system  1,  X  -  exponential^), 

Y  =  time  to  1st  failure  of  system  2,  Y  -  Weibull(  P,  X). 


3Results  obtained  by  solving  equation  (2.1 1)  using  MathCad  3.1,  varying  a  until  the  inequality  was  met. 


15 


Assume:  E[X]  =  -  >  E[Y]  = 

y  px 

Unfortunately  —  <  does  not  imply  that 

Y  P* 

Fx(t)  =  e-tr‘)  >  FY(t)  =  e’u,|B  Vt>C,p>0. 

Therefore,  P(Nx(t)<l)  >  P(NY(t)<l)  cannot  be  shown  using  the  proof  given  for 
Propositions  1  and  2. 

Thus,  the  reader  is  cautioned  not  to  generalize  the  result  obtained  in  the  previous 
analysis.  It  is  unclear  how  the  renewal  process  with  exponentially  distributed  inter  arrival 
times  compares  to  the  renewal  process  with  Weibull  distributed  inter  arrival  times  when 
their  mean  arrival  times  are  equated. 

F.  CONCLUSIONS 

The  preceding  analyses  suggest  that  the  OTA  can  in  fact  reduce  the  amount  of  time 
needed  to  test  a  system  for  mean  time  to  failure  if  the  system  under  test  follows  a  Weibull 
failure  distribution  with  shape  parameter  P  >  1,  and  if  this  shape  parameter  is  known  or 
can  be  estimated.  Unfortunately,  the  shape  parameter  P  is  seldom,  if  ever,  known,  and  it  is 
also  difficult  to  estimate  without  prior  testing,  which  tends  to  make  moot  the  notion  of 
using  the  Weibull  failure  distribution  to  reduce  needed  test  resources.4  So,  from  this 
perspective,  one  can  conclude  that  there  is  probably  little  utility  in  the  OTA  attempting  to 
model  system  reliability  using  the  Weibull  distribution  as  previously  suggested. 

Moreover,  from  the  perspective  of  confidence  levels,  it  remains  unclear  how  the  use 
of  the  test  based  on  the  exponential  distribution,  when  applied  to  a  Weibull  problem, 

4The  Gaver- Jacobs  method  for  "borrowing  information"  from  previous  similar  systems  may  be  useful  in 
this  regard.  Details  of  the  Gaver- Jacobs  method  is  given  in  the  next  section. 
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affects  the  actual  confidence  levels  achieved.  Analysis  suggests  that,  for  purposes  of 
making  inferences  about  the  mean  time  to  failure,  it  might  be  perfectly  safe  for  the  OTA  to 
use  the  exponential  failure  distribution  to  model  a  system  that  would  perhaps  more 
appropriately  be  modeled  using  the  Weibull  distribution,  since,  as  suggested  in  Tables  2 
and  3,  modeling  a  system  which  has  a  Weibull  failure  rate  function  with  the  exponential 
distribution  seems  to  cause  the  OTA  to  overestimate  the  number  of  hours  to  test,  which  in 
turn  seems  to  provide  more  confidence  in  the  result  of  the  test  It  will  be  desirable  to 
investigate  conditions  under  which  this  conservatism  property  is  generally  true. 

Finally,  it  should  be  remembered  that  OTAs  such  as  COTF  tend  to  perform  testing  at 
the  system  level  in  which  there  are  numerous  competing  failure  modes.  It  is  well- 
established  in  the  literature  (e.g.,  NAVORDSYSCOM,  p.  3-25,  1971)  that  as  the  number 
of  competing  failure  modes  increases,  the  occurrence  of  any  one  failure  tends  to  occur  at  a 
constant  rate,  regardless  of  the  individual  underlying  failure  distributions.  So,  despite  the 
warning  concerning  modeling  a  Weibull  system  with  an  exponential  time  between  failure 
distribution,  for  purposes  of  making  inferences  about  the  mean  time  to  failure,  with  larger 
systems  it  is  not  necessarily  inappropriate  to  assume  an  exponential  failure  distribution 
with  a  constant  failure  rate  function. 
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HI.  ON  COMBINING  DT  AND  OT  FAILURE  DATA 


A.  BACKGROUND  AND  PROBLEM  DESCRIPTION 

It  has  recently  been  postulated  that  there  may  be  a  useful  relationship  between 
Developmental  Testing  (DT)  and  Operational  Testing  (OT)  time-between-failure  data  for 
many  systems.  In  fact  it  has  been  noted  anecdotally  by  some  researchers  that  the  OT 
failure  rate  tends  to  be  roughly  four  times  the  DT  failure  rate  for  a  given  system.5  Gaver 
and  Jacobs  (1994)  note  that  if  this  relationship  could  be  quantified,  then  it  could  be  used 
to  anticipate  the  OT  failure  rate  of  a  new,  similar  system  based  on  the  DT  failure  rate  from 
that  system. 

1.  Example  Problem  Setting 

An  OT&E  analyst  has  in  his  possession  historical  time-between-failure  data  for 
both  the  DT  and  OT  phases  of  ten  projects  that  are  similar,  e.g.,  radar  systems.  The 
analyst  believes  a  relationship  exists  between  the  DT  and  OT  failure  rates  such  that 

(i)  =  KAd  (i)  (3.1) 

where  X0  (i)  represents  the  OT  failure  rate  of  system  i,  Xd  (i)  represents  the  DT 
failure  rate  of  system  i,  and  K  represents  a  constant  of  proportionality6  and  where  the 
times  between  failures  are  independent,  having  an  exponential  distribution  with  means 
\/XD  (i)  and  l/Xd  (i) ,  respectively.  The  analyst  further  believes  that  if  he  can  determine  the 

5Interview  between  D.  P.  Gaver,  Professor  of  Operations  Research,  Naval  Postgraduate  School,  and  the 
author,  March  1994. 

6The  author  recognizes  that  this  "constant"  of  proportionality  is  actually  a  random  variable  with  a 
distribution  of  its  own.  Refer  to  section  E  of  this  chapter  for  a  Bayesian  method  of  estimation  of  K. 
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value  of  K,  then  he  can  use  tins  value  on  an  11th  similar  system  to  anticipate  the  OT 
failure  rate  of  System  1 1  given  System  11 's  DT  failure  data. 

B.  THE  GAVER-JACOBS  MODEL 

Suppose  T(j(i;j)  is  a  random  variable  representing  the  time  between  the  (j-l)st  and  the 
jth  time  to  failure  in  DT  of  system  i  where  failures  occur  according  to  a  Poisson  Process 
such  that 


P{Td(i,j)>tJ  =  e^(i)‘‘  ,  td>0. 


(3.2) 


Suppose  further  that  a  similar  relationship  holds  for  T0(i;j),  the  random  variable 
representing  the  time  between  the  (j-l)st  and  the  jth  time  to  failure  in  OT  of  system  i 

P{T0(i,j)>t0}  =  e-x-(i)‘-  ,  to>0.  (3.3) 


Then  Ncj(t(j(i))  is  the  number  of  failures  to  occur  during  DT  for  a  test  of  length  t<}(i) 
fixed  in  advance,  and  has  a  Poisson  distribution,  i.e.: 


P{Nd(td(i))  =  n}  =  e-x<(i)t||(i)  i-d  (l)t<>  (l  - n  =  0,  1, .... 


n: 


(3.4) 


An  equivalent  distribution  holds  for  OT. 

Assume  that  for  each  system  i,  there  are  n^(i)  DT  failures  in  a  test  of  length  tj(i),  and 
n0(i)  OT  failures  in  a  test  of  length  t0(i).  Gaver  and  Jacobs  (p.  3,  1994)  have  devised  a 
maximum  likelihood  estimator  for  the  constant  of  proportionality  K  as  follows: 

The  likelihood  function  is 


L(Xd ,  XD ;  data)  -  e 


(MD) 


n4(i) 


i=l 


nd  (i)! 


n„(i) 


n0(i)! 


(3.5) 
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where  I  denotes  the  total  number  of  similar  systems.  Now  assume  A.0(i)  =  KXd(i);  then, 
up  to  irrelevant  constants,  the  likelihood  of  K  and  the  development  parameters  is 


-^\4(>)t4(i)-IcVx.d(i)t,(i)  i 

L(K,\d(i);data)  =  e  1  ■  ]^[  Xd  (i)“-(i)  K“-(i>X.d  (i)‘ 

i— 1 


n.(i) 


(3.6) 


1=1 


1 

where  nD  =  ^  n0  (i) ,  the  total  number  of  OT  failures  across  all  systems.  Focusing  on  the 

1=1 


parameter  K  above,  take  the  natural  logarithms  of  both  sides  of  equation  (3.6)  to  get 


*(K,  K  (i); data)  =  n0  In  K  -  £  Xd  (i)(td  (i) + Kt0  (i)) 
then  differentiate  with  respect  to  K  to  get 


(3.7) 


i=l 


3K 


n 

K 


(3.8) 


then  finally  set  equation  (3.8)  equal  to  zero  and  solve  for  K  to  get  the  maximum  likelihood 
estimator 


K  =  (3.9) 

Kk 

_  i 

where  Xdt0  =  2^Xd(i)t0(i),  and  ^(i)  is  the  length  of  time  over  which  system  i  is  subjected 

i=l 

to  test  in  OT.  Unfortunately  Xd  (i)  is  unknown,  but  may  be  estimated  as  follows: 

=  (3.10) 

td0) 

The  properties  of  this  estimator  for  K  are  as  yet  unknown  and  will  be  investigated 
shortly.  Meanwhile,  notice  that,  assuming  K  is  well-behaved,  the  analyst  now  has  a 


mechanism  for  relating  past  DT  failure  data  on  a  system  to  future  OT  failure  data  on  the 
same  system  given  historical  DT  and  OT  failure  data  from  similar  systems.  Specifically 


Mi)  =  KXd(i).  (3.11) 

C.  ANALYSIS  OF  THE  GAVER- JACOBS  MODEL 

The  test  of  the  performance  of  the  estimator,  K,  and  the  follow-on  estimate  X0(i), 
was  conducted  using  a  Monte  Carlo  simulation.  The  idea  was  to  simulate  ten  "similar" 
systems  (Systems  1  through  10)  whose  actual  OT  failure  rates  were  four  times  their  actual 
DT  failure  rates.  Each  system  was  tested  for  a  length  of  time,  t,  such  that  the  mean 
number  of  failures  in  (0,  t],  A.t,  was  equal  to  4.0  for  each  system.7  Associated  with  each 
system  was  a  random  number  seed  which  was  used  with  a  Poisson  random  number 
generator,  written  by  the  author  in  Turbo  Pascal  6.0,  the  code  for  which  is  included  in 
Appendix  B.  A  summary  of  the  initial  data  is  given  in  Table  4. 

I.  Estimation  of  K. 

Next,  1000  separate  observations  of  the  number  of  OT  failures  and  DT  failures 
observed  over  the  appropriate  test  time  for  each  system  were  generated  using  the  random 
seed  indicated  in  Table  4.  Then,  K  was  derived  using  equations  (3.9)  and  (3.10)  for  each 
of  the  1000  repetitions  yielding  1000  separate  observations  of  K.  Appropriate  statistics 

-  A 

were  gathered  on  K.  Table  5  summarizes  these  statistics  and  Figure  1  shows  a  histogram 
of  the  observations  of  K. 


7This  number  was  chosen  arbitrarily,  but  was  kept  constant  to  simulate  consistent  choice  of  test  duration. 
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TABLE  4 

INITIAL  DATA  ON  TEN  "SIMILAR"  SYSTEMS  IN  WHICH 


ACT 

[UAL  K  =  4j 

0  AND  ACTUAL  POISSl 

3N  FAILURE  RATE  Xt  =  4.0 

System 

number 

Actual  DT 
failure 

rate 

DT  test 
time 
(hours) 

DT 

Random 

seed 

Actual  OT 
failure 

rate 

OT  test 
time 
(hours) 

OT 

Random 

seed 

I 

0.0002 

20,000 

10 

0.0008 

5,000 

110 

2 

0.0004 

10,000 

20 

0.0016 

2,500 

120 

3 

0.0006 

6,666.67 

30 

0.0024 

1,666.67 

130 

4 

0.0008 

5,000 

40 

0.0032 

1,250 

140 

5 

0.001 

4,000 

50 

0.004 

1,000 

150 

6 

0.002 

2,000 

60 

0.008 

500 

160 

7 

0.004 

1,000 

70 

TiTM 

250 

170 

8 

0.006 

666.67 

80 

166.67 

180 

9 

0.008 

500 

90 

0.032 

125 

190 

10 

0.01 

400 

100 

0.04 

100 

200 

TABLE  5 

SUMMARY  STATISTICS  ON  1000 
OBSERVATIONS  OF  K 


Mean 

4.040 

Std  Error  of  the  Mean 

0.030 

Median 

4.000 

Std  Dev  of  the  simulated  values  of  K 

0.939 

Variance  of  the  simulated  values  of  K 

0.882 

Kurtosis 

0.677 

Skewness 

0.675 

Range 

5.779 

Minimum 

1.913 

Maximum 

7.692 

Sum 

4039.955 

Count 

1000 
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Figure  1  -  Histogram  of  1000  observations  of  estimate  of  K. 

2.  Estimation  of  the  OT  failure  rate. 

While  information  about  the  sampling  distribution  of  K  is  interesting,  the  real 
focus  of  this  study  is  on  how  well  the  analyst  can  predict  the  OT  failure  rate.  Three 
methods  for  estimating  the  OT  failure  rate  are  being  considered  in  this  study: 

1.  Prediction  using  the  relationship  X0(i)  =  KXd(i); 

2.  Direct  observation  of  the  OT  failure  rate  from  OT  data;  and 

3.  Weighted  average  of  the  two  previous  methods. 

In  each  of  the  three  methods,  the  following  are  calculated: 

1 .  The  estimate  for  the  OT  failure  rate; 

2.  The  variance  of  the  estimate,  derived  appropriately  from  each  method; 

3.  The  standard  error  (square  root  of  the  variance)  of  the  estimate;  and 

4.  The  mean-squared  error,  which  is  the  sum  of  the  square  of  the  average  of  the  biases 
and  the  average  of  the  variances  over  all  replications. 


Figure  2,  below,  describes  one  replication  of  the  OT  failure  rate  estimation 


simulation  model. 


Figure  2  •  One  replication  of  the  OT  failure  rate  estimation  simulation  model. 


a.  Prediction  using  the  estimator  for  K. 

The  analyst's  prediction  of  the  future  OT  failure  rate  for  an  11th  system, 
using  the  estimator  for  K,  would  be 

X.(11)-KM11).  (3.12) 

The  associated  variance  would  be  computed  as  follows: 

Var[X0(ll)J  =  Var[KXd  (11)] 

=  E[(iad(ll)-E[K]E[\d(ll)])2] 

=  E[K2Xd  (1 1)  2  -  2KXd  (1 1)  E[K]E[Xd  (11)]+  E[K]2£[Xd  (1 1)  ]2] 

=  E[K2]E[Xd  (1 1) 2]  -  E[K]2£[X.d  (1 1) ]2] 

=  (Var[K]  +  E[K]2  )(Var[Xd(ll)]  +E[Xd(ll)]2  )  -  E[K]2E[Xd(ll)]2 
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simplifying  gives: 

Varied  1)]  =  Var(K]Var[A.d(l Dl  +  Var[K]E[Xd (1 1)]2  )  +  Var[Xd(ll)]E[K]2  .  (3.13) 

A  A 

Expressions  are  needed  for  Var[K]  and  Var[Xd(ll)].  An  approximation  to 

a 

Var[K]  can  be  derived  through  maximum  likelihood  estimation  using  observed  Fisher 
information  as  follows  (Cox  and  Hinkley,  p.297,  1974): 


Var[K]  = 


(3.14) 


3KMk=K 


Taking  the  second  derivative  with  respect  to  K.  of  equation  (3.7)  gives 


d2e  _  n0 


9K2 


K 


2  ' 


(3.15) 


Thus 


K2 

Var[K]  =  — 
n 


n„ 


_  V 


n 


o  __  _ 


(W1 


St 


nd(i) 
(i) 


(3.16) 


t.(i) 


Var[Xd  (11)]  can  be  derived  as  follows  by  noting  that  Xd  (i)  is  simply  a  Poisson  random 
variable,  Ncj(i),  multiplied  by  a  constant,  * 


td(i) 


Var[A.d(i)]  =  Var 


(td(i)) 


f  Nd(i) 

[tdO)J  (td(i)): 


1  r\d(i)td(i)«-^l 


Var[Nd(i)) 


(3.17) 


(td(i))2 
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The  standard  error  of  Xd  (1 1)  is  then  simply  the  square  root  of  Var[Xd  (1 1)  ]. 

An  example  prediction  can  now  be  performed.  Suppose,  using  data  from  the 

a 

first  simulation  run,  equation  (3.9)  gives  K  =  4.293  and  equation  (3.10)  gives  Xd(ll)  = 
0.0003.  Then 

ia  (1 1)  =  KXd  (1 1)  =  (4.293X0.0003)  =  0.00129. 


Note  that  per  the  original  assumption,  K  =  4.0,  so  we  would  expect  a  system  whose  DT 
failure  rate  is  0.0003  to  have  an  OT  failure  rate  of  0.0012.  So,  at  first  glance  this  appears 
to  be  a  relatively  good  estimate  of  the  true  situation. 

The  variance  and  standard  error  of  X,0(ll)  is  solved  using  equation  (3.13)  as 

follows: 


E[K]  =  K=  4.293 
E[A,d(ll)]  =  A.d  (1 1)  =0.0003 


Var[K]  =  - — 

(  10 

I 

V  i=l 

Var[Xd(ll)]  = 


n„ 


-  no  - 


44 


nd(i) 

tdO) 


t0(i) 


nd  (1 1) 


6 


(td(ll))2  (20.000)2 


(10.25) 


=  1.5x10' 


0.4188 


finally 


Var[X.0(ll)]  =  Var[KXd(ll)]  =  (0.4188)(1.5xl0'8)  +  (L5xl0'8)(4.293)2 
+  (0.4 188)(0.0003)2  =  32xl0'8. 

And  the  standard  error  of  X0(ll)  is  V32  x  10-8  =  5.66  x  10"^.  Table  6  summarizes  the 
results  of  this  prediction. 
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TABLE  6 

SUMMARY  OF  RESULTS  OF  PREDICTION  OF  OT  FAILURE  RATE 


ON  SYSTEM  1 1  USING  ESTIMATE  O 

F  K  BASED  ON  FIRST  REPLICATION 

Actual  Failure  Rate  X.0  (1 1) 

0.0012 

Estimated  Failure  Rate  (1 1) 

0.00129 

Variance[X0(ll)] 

32xl0-8 

Standard  Error[A,n  (1 1)  ] 

5.66  x  10"4 

X.0(ll)  ±  1  Standard  Error 

( 0.0007,  0.0019 ) 

Xn  (1 1)  ±2  Standard  Errors 

( 0.0002,  0.0024 ) 

b.  Direct  observation  from  OT  failure  data. 

A. 

Of  course,  the  analyst  could,  if  he  desired,  come  up  with  a  value  for  XG(1 1) 
by  observing  the  empirical  OT  failure  rate  as  usual.  In  this  situation,  the  analyst  would 
note  the  number  of  failures  which  occur,  n^ll),  during  the  time  on  test,  t0(ll),  and 
derive  the  OT  failure  rate  as  follows: 


x;ou 


n„(ll) 

t0(H) ' 


(3.18) 


For  example,  suppose  that  System  1 1  is  operationally  tested  for  a  period  of 
5,000  hours  and  experiences  6  failures  during  that  time.  Then  the  OT  failure  rate 
observed  through  operational  testing  would  be 


*:ai) 


=  — - —  =  0.0012. 

5000 


Using  the  same  reasoning  as  in  equation  (3.17),  the  variance  of  the  observed 
value  is  calculated  as 
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Var[X'(l  1)]  =  Var 


'no(H)‘ 

r  i  \ 

_t0(ll) 

U(11)J 

Var(N0(l  1)) 


'  1  V 


t0(ll) 

6 


X0(ll)t0(ll)  = 
=  24x10"* . 


tQ(ll) 


n0(l  1) 


v-o 


(5000)2 

And  the  standard  error  of  X'(ll)  is  V24x  10~8  =  4.9  x  10*4.  Table  7  summarizes  the 
results  of  these  calculations. 


TABLE  7 

SUMMARY  OF  RESISTS  OF  DIRECT  OBSERVATION  OF  OT  FAILURE  RATE 
ON  SYSTEM  1 1  FROM  OT  FAILURE  DATA 


Actual  Failuic  Rate  XQ  (1 1) 

0.0012 

Estimated  Failure  Rate  X'  (1 1) 

0.0012 

Variance[X'(ll)] 

24  x  10’8 

Standard  Eiror[X' (1 1)3 

4.9  xlO*4 

X'(ll)  ±  1  Standard  Error 

(0.0007, 0.0017) 

X'  (11)  ±2  Standard  Errors 

(0.0002, 0.0022) 

Unfortunately,  this  method  would  not  save  the  analyst  any  time  since  he 
would  have  to  run  the  entire  operational  test  procedure  to  derive  a  meaningful  number. 
Moreover,  he  would  not  be  taking  advantage  of  the  available  DT  failure  data, 
c.  Weighted  average  of  predicted  value  and  direct  observation. 

Suppose  now  that  the  analyst  runs  the  entire  operational  test  procedure  as 
above,  but  also  desires  to  use  the  available  DT  failure  data  in  an  effort  to  improve  his 
estimate  for  the  OT  failure  rate.  Assuming  the  analyst  has  the  data  from  the  10  other 
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previous  systems  with  which  to  derive  the  prediction,  X0(ll)  =  KXd(ll),  he  can  then 
calculate  a  weighted  average  of  the  two  estimates. 

An  estimator  for  the  mean  of  a  linear  combination  of  two  random  variables, 
X  and  Y,  with  common  mean  m  and  different  variances,  a*  and  Gy  can  be  constructed  as 
follows: 

Let  m  =  aX  +  bY ,  0<a,b£i,  a  +  b=l. 

Then  E[m]  =  am  +  bm  =  m(a  +  b)  =  m(l)  =  m, 

*(a)  =  E[(m-  m)2]  =  E[(aX  +  (l-a)Y  -  m)2] 

=  E[(a(X  -  m)  +  (1  -  a)(Y  -  m))2] 

=  E[a2(X-m)2  +  2a(l-a)(X-m)(Y-m)  +  (I-a)2(Y-m)2] 

=  a2Var[X]  +  (1  -  a)2Var[Y]- 

Differentiating  with  respect  to  a  and  setting  equal  to  zero  gives 

~ -  =  2aVar[X]  -  2(1  -  a)Var[Y]  =  0 
da 

which  yields 


a[Var[X]  +  Var[Y]J  =  Var[Y] 


or 


a  _  Var(Y] 

Var[X]  + VarlY]  ' 


Similarly 


Var[X]  Xt 

Var[X]  +  Var[Y]  /al  +  /a>' 
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Substitution  then  yields 


A 


(3.19) 


Therefore  an  estimator  for  the  OT  failure  rate,  XQ(i),  can  be  constructed  as  follows: 


K(i) 


Mi>  = 


KXdU)  _ 

Var[idd(i)]'rVar[A.:(i)] 


(3.20) 


Var[KXd(i)j  Var[X;(i)] 


A  A  * 

where  KAd  (i)  represents  the  predicted  value  of  the  OT  failure  rate  and  A  '(i)  denotes  the 

empirical  OT  failure  rate  observed  during  operational  testing. 

For  example,  the  values  obtained  previously  for  KAd(ll),  the  predicted 
value,  and  A.' (11),  the  observed  value,  were  0.00129  and  0.0012  respectively.  It  was 
shown  earlier  that  the  estimated  variance  for  the  predicted  value  was  Var[KAd  (1 1)]  =  32 
x  10*8,  and  the  variance  of  the  observed  value  was  Var[X.'(ll)]  =  24  x  10*8.  Inserting 

these  numbers  in  equation  (3.20)  gives 

12.9  x  IQ-4  12x10"* 

i0  (1 1)  =  32x^10:! - ?4xl£l  - 12. 39  x  10"*. 

32x10"®  +  24x10"® 

The  variance  associated  with  this  estimate  can  be  constructed  as  follows: 

Let  Var[m]  =  Var[aX  +  bY]  =  a2Var[X]  +  b2Var[Y]. 
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Substitution  yields 


oi  = 


(/J2<*x+(/J2° 


(3.21) 


Therefore,  the  estimated  variance  of  X.0(l  1)  is  computed  as  follows: 


Vara  (11)]  =  • 


1  1 
f - 


1  1 
-+- 


=  13.7x10  . 


Var[KA.d(ll)]  Var[X^(ll)]  32x10"*  24x10' 


And  the  standard  error  of  Xn(ll)  is  Vl3.7  x  10"*  =  3.7  x  10'4.  Table  8  summarizes  the 
results  of  these  calculations. 


TABLE  8 

SUMMARY  OF  RESULTS  OF  WEIGHTED  AVERAGE  OF 


PREDICTED  VALUE  AND  DIRECT  01 

BSERVATION  OF  OT  FAILURE  RATE 

Actual  Failure  Rate  Xo(ll) 

0.0012 

Estimated  Failure  Rate  X0  (1 1) 

0.00124 

Variance[X0(ll)] 

13.7  x  10-8 

Standard  Error[X,c  (1 1)  ] 

3.7  x  10-4 

i0  (11)  ±  1  Standard  Error 

(0.0009,  0.0016) 

(11)  ±2  Standard  Errors 

( 0.0005, 0.0020 ) 

3.  Comparison  of  the  three  methods  for  estimating  the  OT  failure  rate. 

Table  9  summarizes  the  results  obtained  for  each  of  the  three  methods  of  deriving 
the  OT  failure  rate  for  one  replication  of  the  simulation.  They  are  arguably  quite  similar. 
So  to  answer  the  question  of  which  method  is  best,  further  analysis  is  required. 
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TABLE  9 

COMPARISON  OF  THE  THREE  METHODS  FOR  OBTAINING 
THE  OT  FAILURE  RATE 


Statistic 

Prediction  Method 

Direct  Observation 
Method 

Weighted  Average 
Method 

Actual  Value 

0.0012 

0.0012 

0.0012 

Estimated  Value 

0.00129 

0.0012 

0.00 ’24 

Variance 

32  x  10-8 

24  x  10'8 

13.7  x  10-8 

Standard  Error 

5.66  x  10-4 

4.9  x  10'4 

3.7  x  10-4 

Estimate  ±  1  S.E. 

( 0.0007, 0.0019 ) 

( 0.0007,  0.0017  ) 

( 0.0009, 0.0016 ) 

Estimate  ±  2  S.E.s 

( 0.0002,  0.0024 ) 

( 0.0002, 0.0022 ) 

( 0.0005, 0.0020 ) 

Since  the  results  of  the  three  methods  discussed  previously  were  so  similar,  a 
look  at  the  asymptotic  behavior  of  the  three  methods  is  in  order.  To  facilitate  this 
analysis,  four  new  "similar"  systems  were  created,  systems  11,  12,  13  and  14,  each  with 
failure  rates  and  test  times  of  similar  construction  to  the  previous  10  systems,  i.e., 
X0(i)  =  KXd(i),  where  K  is  equal  to  4.0  and  Xt  is  also  equal  to  4.0.  Each  system  also  had 
its  own  random  seeds  for  both  DT  and  OT  failures.  Table  10  summarizes  the  four  new 
systems. 


TABLE  10 

INITIAL  DATA  ON  FOUR  NEW  "SIMILAR"  SYSTEMS  IN  WHICH 
ACTUAL  K  =  4.0  AND  ACTUAL  POISSON  FAILURE  RATES  Xt  -  4.0 


System 

number 

Actual  DT 
failure 

rate 

DT  test 
time 
(hours) 

DT 

random 

seed 

Actual  OT 
failure 
rate 

OT  test 
time 
(hours) 

OT 

random 

seed 

11 

0.0003 

13,333.33 

12 

0.0012 

3,333.33 

112 

12 

0.0005 

8,000 

22 

0.002 

2,000 

122 

13 

0.0007 

5,714.29 

32 

0.0028 

1,428.57 

132  ! 

14 

0.0009 

4,444.44 

42 

0.0036 

1,111.11 

142 
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One  thousand  repetitions  of  a  simulated  OT  and  DT  were  run  on  each  of  the  four 
new  systems,  and  then  the  exact  same  analysis  as  that  done  previously  in  section  III.C.2 
was  performed  for  each  replication.  Grand  estimates  appear  in  Tables  11  through  14  in 
the  column  entitled  "Average  value  (n  =  1000)"  and  were  produced  by  averaging  the 
estimates  of  the  following  from  all  1000  simulated  observations,  or  sub-estimates: 

1 .  The  estimated  OT  failure  rate  derived  from  prediction,  direct  observation  and 
weighted  average; 

2.  The  estimated  variance  of  the  estimated  OT  failure  rate  derived  from  prediction, 
direct  observation  and  weighted  average;  and 

3.  The  estimated  standard  error  of  the  estimated  OT  failure  rate  derived  from 
prediction,  direct  observation  and  weighted  average. 

Then,  in  order  to  get  an  idea  of  the  sample  distribution  of  each  of  the  above 
super-estimates,  the  sample  variance,  sample  standard  deviation  and  coefficient  of 
variation8  were  computed  for  the  1000  sub-estimates  for  each  system.  Finally,  the 
estimated  mean-squared  error,  which  is  the  sum  of  the  square  of  the  average  of  the  biases 
with  the  average  of  the  variances  over  all  replications,  i.e., 

MSF  =  (Average  of  the  biases)2  +  (Average  of  the  variances), 
was  calculated  for  each  prediction  method  and  placed  in  the  last  column  of  each 
table.  Figure  3  describes  the  situation. 


^^The  coefficient  of  variation  is  simply  the  sample  standard  deviation  divided  by  the  sample  average,  and 
serves  as  a  measure  of  the  spread  of  the  sample. 
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Figure  3  -  Methodology  for  comparing  different  ways  of  estimating  OT  failure  rate. 


I 
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a.  New  System  11. 


The  statistics  on  new  System  1 1  are  summarized  in  Table  1 1  and  a  graphical 
summary  is  shown  in  Figure  4.  Recall  that  actual  (1 1)  =  0.0012. 


TABLE  11 

NEW  SYSTEM  11  -  ESTIMATED  OT  FAILURE  RATES  AND  OTHER 
RELEVANT  STATISTICS  (ACTUAL  Xn(ll)  =  0.0012) 


Prediction  Direct  Observation  Weighted  Average 


Figure  4  -  Graphical  summary  of  System  11  estimates. 


b.  New  System  12. 


The  statistics  on  new  System  12  are  summarized  in  Table  12  and  a  graphical 
summarv  is  shown  in  Figure  5.  Recall  that  actual  (12)  =  0.002. 


TABLE  12 

NEW  SYSTEM  12  -  ESTIMATED  OT  FAILURE  RATES  AND  OTHER 


Sub¬ 

estimates 

HI 

Coeffof 

Variation 

Mean 

Squared 

Error 

Prediction 

0.0020 

KfiSEH 

0.0011 

55% 

ISH 

msm 

70% 

1.23X10"6 

EUtBn 

0.0010 

mmam 

0.00038 

38% 

Direct 

0.0020 

jHHHM 

0.00098 

49% 

Observation 

mi 

WSZM 

49% 

9.99X10’7 

0.000% 

7.42x10"^ 

0.00027 

28% 

Weighted 

Mjrnm 

0.0017 

0.00077 

45% 

Average 

V[i.(12)l 

man 

49% 

5.49xl0'7 

0.00065 

3.71x10'^ 

0.00019 

29% 

I  Actual  Failure  Rata 
0.002 


New  System  12 

I  Mean  Estimated  Failure  Rats  BMsan  Estimated  Std  Error 


Prediction 


Direct  Observation  Weighted  Average 


Figure  5  -  Graphical  summary  of  System  12  estimates. 
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c.  New  System  13. 


The  statistics  on  new  System  13  are  summarized  in  Table  13  and  a  graphical 
summary  is  shown  in  Figure  6.  Recall  that  actual  X0  (13)  =  0.0028. 


TABLE  13 

NEW  SYSTEM  13  -  ESTIMATED  OT  FAILURE  RATES  AND  OTHER 
RELEVANT  STATISTICS  (ACTUAL  Xn(13)  =  0.0028) 


d.  New  System  14. 

The  statistics  on  new  System  14  are  summarized  in  Table  14  and  a  graphical 
summary  is  shown  in  Figure  7.  Recall  that  actual  X0(14)  =  0.0036. 


TABLE  14 

NEW  SYSTEM  14  -  ESTIMATED  OT  FAILURE  RATES  AND  OTHER 


Prediction  Direct  Observation  Weighted  Average 


Figure  7  -  Graphical  summary  of  System  14  estimates. 
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D.  RESULTS  OF  THE  GAVER-JACOBS  MODEL  ANALYSIS 

Simulation  results  indicate  that  the  predicted  value  method  and  the  observed  value 
method  are  unbiased  and  provide  good  estimates  of  the  OT  failure  rate.  The  weighted 
average  method,  although  it  gives  a  more  precise  estimate  because  of  its  smaller  mean 
squared  error,  is  biased  in  that  it  tends  to  underestimate  the  actual  failure  rate.  The 
observed  value  method  gives  slightly  better  precision  than  the  predicted  value  method, 
since  the  mean  squared  error  obtained  with  the  observed  value  method  tends  to  be  slightly 
less  than  that  associated  with  the  predicted  value  method.  Notice,  however,  that  except 
for  the  very  small  price  paid  in  terms  of  precision,  the  estimate  obtained  using  the  Gaver- 
Jacobs  predicted  value  method  seems  to  estimate  the  OT  failure  rate  nearly  as  well  as 
direct  observation  of  the  OT  failure  rate.  This  is  significant  since  the  Gaver-Jacobs 
estimate  can  be  used  prior  to  the  start  of  operational  testing  (provided  the  constant-K 
model  is  nearly  correct). 

The  bottom  line  is  that  if  the  relationship  between  the  DT  failure  rate  and  the  OT 
failure  rate  proposed  earlier  does  in  fact  exist,  then  the  Gaver-Jar  'bs  method  of  estimating 

A  A  A 

the  OT  failure  rate  of  a  new  system,  X0(i)  =  KXd(i),  given  that  system's  DT  failure  rate  as 
well  as  DT  and  OT  data  on  previous  similar  systems,  seems  to  provide  a  reasonable 
estimate  of  the  OT  failure  rate  of  the  new  system  which  performs  nearly  as  well  as  direct 
observation  of  the  OT  failure  rate  through  operational  testing.  Moreover,  the  Gaver- 
Jacobs  estimate  of  the  OT  failure  rate  can  be  computed  prior  to  the  start  of  operational 
testing. 
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E.  A  HIERARCHICAL  BAYES  MODEL  FOR  ESTIMATION  OF  K 

It  is  idealistic,  in  the  previous  model,  to  assume  that  all  systems  have  exactly  the  same 
value  for  K.9  Rather,  assume  that  each  system  i  has  its  own  K  value,  Kj,  and  that  it  is 
selected  independently  from  a  prior  distribution.  The  gamma  distribution,  for  example, 
turns  out  to  be  conjugate,  so  integrals  can  be  explicitly  computed.  Its  density  function  is: 

qK(K)  =  e-'*l£^-p.  (3.22) 


From  equations  (3.5)  and  (3.6)  it  can  be  seen  that  the  likelihood  of  Kj,  up  to  irrelevant 
constants,  is 

L,  (Kj;  Xd  (i);i  -data)  =  e'IE.^(i)t‘<i)K."'(i>[e-i'(i),-(i)(Xd  (i))n,,(i)+n-(i)].  (3.23) 


To  obtain  the  marginal  posterior  density  of  Kj,  equation  (3.22)  is  re-written,  up  to 
irrelevant  constants,  as  follows: 


qK  (K,;l- data)  =  - (p + *a ©t. (i)).  (3.24) 

l(0+no(i)) 


If  the  parameters  p,  8  and  Xd  (i)  were  known,  and  the  data  available,  i.e.,  n0(i)  and 
t„(i),  then  one  could  compute  the  conditional  expected  value  of  Kj  and  use  this  as  a  point 
estimate.  However  there  are  other,  more  interesting  ways  to  use  the  posterior.  First,  it  is 


9The  following  explanation  for  the  hierarchical  Bayes  model  was  developed  by  D.  P.  Gaver,  Professor  of 
Operations  Research,  Naval  Postgraduate  School,  in  handwritten  notes  to  the  author,  July  1994. 


necessary  to  evaluate  the  various  parameters  p,  8  and  Xd(i).  To  do  this,  the  "empirical 
Bayes"  method  is  used,  meaning  that  the  unknown  Kj  term  is  integrated  out  of  equation 
(3.23)  to  obtain 


*, *.<#.«>  K..<i>i^.(pKi)i-,p  dKi 

1  (o) 


(3-25) 


The  integral  in  the  above  expression  can  be  reduced  to 


_pL 

r(8) 


(n.(i)+6)-l  -0^(i)t.(i)+f>)K, 


0 


dK. 


(3.26) 


Noting  that  the  cumulative  distribution  function  (CDF)  of  a  gamma  distribution  with 
parameters  (nc(i)  +  8)  and  (Xd(i)t0(i)  +  p)  has  a  known  total  area  of  1.0,  the  following 
identity  can  be  stated: 


(Xd(i)t0(i)  +  p) 
r(n0(i)  +  5) 


n.(i)+S  «, 


'  J  K  (i)t.  (i)*^p)X,  ^  j 


(3.27) 


Rearranging  equation  (3.26)  to  take  advantage  of  (3.27),  equation  (3.26)  can  be  reduced 


to 


r(n0(i)  +  S)p8 


r(5)(Xd(i)t0(i)  +  p) 


(3.28) 


Therefore,  equation  (3.25)  can  be  re-written  as 


L,  (Xd  (i);  p;  8;  i  -  data)  =  [e'x'‘(,),d  (i)  (Xd  (i)) 


n.(i)+na(i) 


r(n  (i)  +  S)p‘ 


r(8XX,(i)t.(i)+p) 


n,(i)+8 


.(3.29) 


It  is  generally  easier  to  work  with  the  natural  log  of  the  likelihood  function 


£i  =  ln(L, )  =  -Xd  (i)td  (i) + (nd  (i)  +  nc  ri))ln(Xd  (i))  +  ln(r(nG  (i)  +  8)) 
+  5ln(p)  -  ln(r(8))  -  (n0  (i)  +  8)ln(Xd  (i)t0  (i)  +  p) 


(3.30) 


The  partial  derivative  of  equation  (3.30)  with  respect  to  Xd  (i)  is 


dii 

axda) 


nd(i)-Xd(i)td(i)  t 

Mi) 


n0(i)--Xd(i)t0(i) 
P _ 

Mi)(Xd(i)+p) 


(3.31) 


Setting  equation  (3.31)  equal  to  zero  and  solving  for  Xd  (i)  gives  the  maximum  likelihood 
estimator  for  Xd  (i).  Unfortunately,  however,  the  parameters  8  and  p  are  unknown. 

One  way  to  get  at  the  parameters  8  and  p  is  to  obtain  an  initial  estimate  for  Xd  (i) 
using  the  formula 


(i)  =  ~“»  (3-32) 

td  (i) 

and  then  derive  maximum  likelihood  estimators  for  8  and  p.  One  can  then  use  these  initial 
estimates  for  8  and  p  in  equation  (3.31)  to  obtain  an  adjusted  value  for  Xd(i)  for  each 

A  A 

system.  It  may  be  possible  to  continue  these  iterations  until  8 ,  p ,  and  Xd  (i)  all  stabilize. 


Then,  using  8  and  p ,  as  the  parameters  for  the  gamma  prior  distribution,  it  is  possible  to 
estimate  the  mean  Kj  factor  and  its  variance  using  the  relationships 


and 


(3.33) 


(3.34) 


Assuming  data  exists  on  more  than  one  "similar"  system,  say  I  of  them,  then  the  log- 
likelihood  function  becomes,  up  to  irrelevant  constants 


i=l 


Sln(p)+ln 


r(n0(i)+S) 

r(8) 


-  (n0  (i)  +  8)  ln(kd  (i)t„  (i)  +  p) 


.  (3.35) 


The  ratio  of  gamma  functions  in  equation  (3.35)  can  be  simplified  by  expanding  the 
numerator  and  then  canceling  the  denominator 

r(n„(i)+8)  (n.  (i ) + 5  -  l)(n.  (i) + 5  -  2)-  -  -  (n0  (i ) + S  -  n0  (i  ))r(S) 

r(s)  r(8) 

(1,  no(i)  =  0  (3-36) 

|(n0 (i) + 8 - l)(nc (i) + 8 - 2)- •  -(nc (i) + 8 - nc (i)),  n0(i)>l. 


Then  the  natural  log  of  equation  (3.36)  becomes 

0,  nQ(i)  =  0 

n,(i)-l 

2^1n(8+ j), 

.  j=° 


f‘(8)  =  ln 


r(n0(i)  +  8) 

r(S) 


n0(i)>l. 


(3.37) 
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Rewriting  equation  (3.35)  gives 


^  =  £[8ln(p)  +  f*  (8)-{n0(i)+8)ln(\d(i)t0(i)  +  p)] 

i=l 


(3.38) 


which  can  be  further  simplified  to 

1  -  £ [8  ln(p)  - (n0 (i)  +  8) ln(Xd (i)t0  (i)  +  p)]  +  £  f '(5) . 
1=1  1*1 


(3.39) 


It  is  now  a  simple  matter  to  take  the  partial  derivatives  of  equation  (3.39)  with  respect  to 
8  and  p  which,  when  set  equal  to  zero,  will  enable  one  to  solve  for  maximum  likelihood 
estimators  for  the  two  parameters.  Specifically 


where 


and 


di 

35 


Iln(p)-j>(Xa(i)Mi)+p)+X^f(8) 

i=l  "O 
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f*(S)  = 


fO,  no(i)  =  0 

It- .  n0mi 

M  8  +  j 


dp 


j  -KdK(i)-n0(i) 

=  y  p _ 

tr  (Mi)t0(i)+p) 


(3-40) 


(3.41) 


(3-42) 


Finally,  the  preceding  two  equations,  (3.40)  and  (3.42)  along  with  (3.31)  give  1  +  2 
equations  with  1  +  2  unknowns  which  can  be  solved  to  obtain  Xd  (i) ,  8  and  p . 

So  far,  the  author  has  not  been  able  to  solve  the  system  of  1+2  equations  and  1+2 
unknowns.  Solution  with  1  =  10  similar  systems  was  attempted  using  the  non-linear  solver 
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included  with  Microsoft  Excel  4.0  with  solver  options  set  as  follows:  i.^wion  search 
method,  using  forward  derivatives,  with  tangential  estimates.  Moreover,  attempts  to  use 
the  same  methodology  with  only  two  "similar"  systems;  i.e.,  four  v-quauons  and  four 
unknowns,  have  also  failed.  The  author  believes  that  equations  (3.31),  (3.40)  and  (3.42), 
or  at  least  some  combination  of  the  equations,  are  ill-conditioned  in  that  they  fail  to 
converge  at  a  common  point  in  (1+2)  space.  An  attempt  was  made  to  re-parameterize  the 
equations,  substituting  the  variables  Xd(i),  p  and  5  with  exponentials,  however  this 
attempt  also  failed  to  convene. 

While  it  appears  the  hierarchical  Bayes  methodology  has  promise,  as  well  as  a  certain 
intuitive  appeal,  further  research  is  called  for  in  order  to  better  condition  the  set  of 
equations  needed  to  use  the  model.  Other  issues  to  consider  are  whether  or  not  the 
gamma  prior  is  appropriate  for  this  model  and  if  so,  how  the  set  of  equations  needed  to 
use  the  model  can  be  simplified  or  re-parameterized  so  that  reliable  model  parameter 
estimates  can  be  obtained. 
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IV.  SUMMARY,  CONCLUSIONS  AND  RECOMMENDATIONS 

A.  SUMMARY 

The  purpose  of  this  thesis  is  to  ( 1 )  investigate  the  behavior  of  a  statistical  procedure 
used  to  determine  the  length  of  time  needed  to  test  a  continuously  operating  system  when 
the  system  under  test  is  assumed  to  have  an  exponential  time  to  failure  distribution,  but 
when  in  fact  the  system  does  not  have  an  exponential  time  to  failure  distribution  and  rather 
may  have  a  time  to  failure  distribution  that  more  closely  resembles  a  Weibull  distribution; 
(2)  to  examine  the  behavior  of  the  Gaver-Jacobs  estimator  for  the  relationship  between 
DT  and  OT  failure  rates  using  a  Monte  Carlo  simulation;  and  (3)  to  introduce  a 
hierarchical  Bayes  method  for  estimation  of  the  relationship  between  the  DT  and  OT 
failure  rates.  In  the  first  pan  of  the  thesis,  the  parameters  of  the  procedure  are  chosen 
based  or.  the  assumption  of  an  exponential  time  to  failure  distribution.  The  behavior  of  the 
procedure  is  then  investigated  under  the  assumption  that  the  time  between  failures  are 
independent,  having  a  Weibull  distribution  with  the  same  mean  time  to  failure.  The 
intention  is  to  test  whether  or  not  use  of  the  Weibull  distribution  will  result  in  fewer  test 
hours  to  verily  that  a  system  under  test  has  a  mean  time  to  failure  which  meets  or  exceeds 
a  minimum  acceptable  value  at  a  given  level  of  confidence,  and  also  to  test  whether  or  not 
the  O'I'As  are  misstating  the  confidence  they  report  in  their  operational  testing  results 
when  the  exponential  distribution  is  used  to  model  a  system  whose  true  failure  rate 
function  is  that  of  a  Weibull  distribution.  In  the  second  part  of  the  thesis,  the  Gaver- 
Jacobs  estimator  is  examined  in  terms  of  the  consistency  and  variability  of  the  OT  failure 
rate  derived  using  that  estimator  compared  to  other  methods  of  establishing  the  OT  failure 
rate.  Finally,  in  the  third  part  of  the  thesis,  a  hierarchical  Bayes  model  is  proposed  for  the 
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ultimate  purpose  of  being  able  to  describe  the  mean  and  variance  associated  with  a  random 
relationship  between  DT  and  OT  failure  rates. 

B.  CONCLUSIONS 

1.  On  comparison  of  the  exponential  and  YVeibull  failure  models. 

The  first  analysis  suggests  that  operational  testing  agencies  can,  under  certain 
circumstances,  reduce  the  amount  of  time  needed  to  demonstrate  that  the  mean  time  to 
failure  (MTTF)  of  a  continuously  operating  system  is  at  least  as  large  as  the  minimum 
acceptable  value  (MAV)  if  the  system's  true  underlying  failure  distribution  is  Weibull  with 
shape  parameter  |3  >  1 .  However,  the  shape  parameter  is  seldom,  if  ever,  known  or  even 
estimable  without  expending  valuable  test  resources,  so  it  is  unlikely  that  an  OTA  could 
actually  save  any  time  or  resources  with  this  methodology.  Moreover,  from  the 
perspective  of  confidence  levels,  the  analysis  suggests  that  confidence  in  one's  result  is  not 
sacrificed  when  the  exponential  distribution  is  used  to  model  a  system  which  has  a  non¬ 
constant  failure  rate,  since  it  appears  that  use  of  the  exponential  distribution  on  a  truly 
Weibull  system  seems  to  cause  the  OTA  to  overestimate  the  number  of  hours  needed  to 
test  the  MAV  for  the  MTTF  to  a  given  confidence  level,  and  also  to  underestimate  the 
confidence  it  reports  in  its  results  under  the  incorrect  assumption  of  an  exponential 
distribution.  However,  for  reasons  mentioned  previously,  these  results  may  not  be  true  in 
general.  Therefore,  it  is  concluded  that  there  is  currently  little  to  be  gained  by  attempting 
to  model  the  time  between  system  failures  with  the  Weibull  distribution  as  previously 
suggested,  particularly  if  very  few  observations  of  time  to  failure  are  available. 

2.  On  using  the  Gaver-Jacobs  method  to  combine  DT  and  OT  failure  data. 

The  second  analysis  indicates  that  the  Gaver-Jacobs  method  for  estimating  a 

system's  OT  failure  rate  given  that  system's  DT  failure  rate  as  well  as  failure  data  on 
previous  similar  systems  produces  a  reasonable  estimate  of  the  actual  OT  failure  rate 
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(when  the  OT  failure  rate  behaves  as  expected).  In  other  words,  if  the  relationship 
between  the  DT  and  OT  failure  rates  for  a  certain  system  is  such  that  it  can  in  fact  be 
estimated  from  previous  similar  systems,  then  the  Gaver-Jacobs  method  gives  a  good 
estimate  of  the  OT  failure  rate  given  that  system's  DT  failure  rate.  In  fact,  the  estimate 
provided  by  the  Gaver-Jacobs  method  appears  to  be  nearly  as  good  as  direct  observation 
of  the  OT  failure  rate  obtained  through  actual  testing.  Thus,  if  it  can  be  shown  that 
inferences  about  the  failure  rate  of  a  system  can  in  fact  be  made  on  the  basis  of  previous 
similar  systems,  then  the  Gaver-Jacobs  method  of  predicting  a  new  system's  OT  failure 
rate  from  previous  DT  failure  data  has  some  potential  utility.  In  particular,  this  OT  failure 
rate  prediction  could  be  used  in  ways  which  would  aid  decisions  as  to  when  to  begin  OT 
as  well  as  to  augment  follow-on  OT  in  an  effort  to  reduce  needed  testing  resources. 

3.  On  the  hierarchical  Bayes  Model  for  estimating  K. 

The  third  analysis  is  incomplete.  It  was  hoped  that  the  hierarchical  Bayes  method 
given  would  provide  a  means  by  which  to  estimate  K,  the  relationship  between  the  DT  and 
OT  failure  rates,  when  K  is  assumed  to  be  a  random  variable  with  some  prior  distribution. 
Regretfully,  the  equations  derived  for  the  maximum  likelihood  estimators  for  the 
parameters  of  the  gamma  prior  have  proved  to  be  difficult  to  solve  satisfactorily.  The 
author  was  unable  to  find  a  meaningful  solution  to  the  set  of  1+2  equations  and  1+2 
unknowns,  where  I  represents  the  number  of  "similar"  systems  for  which  previous  DT  and 
OT  data  exist.  If  such  a  solution  were  available  it  could  be  used  to  predict  future  failure 
rates  by  means  of  a  Bayesian  calculation.  Investigation  of  such  a  procedure  is  left  for  a 
later  time. 

C.  RECOMMENDATIONS  FOR  FURTHER  STUDY 

Concerning  the  notion  of  modeling  the  failure  distribution  of  a  continuous  system 
using  the  Weibull  distribution,  OTAs  could  clearly  benefit  from  being  able  to  use  this 
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approach  when  appropriate.  If  the  time  between  failures  for  a  system  has  a  Weibull 
distribution,  then  proper  choice  of  that  system's  Weibull  shape  parameter  would  probably 
enable  the  OTA  to  reduce  the  required  test  resources  while  maintaining  a  specified  level  of 
confidence  in  the  result  of  their  test.  Unfortunately,  there  is  very  little  data  in  existence 
which  can  assist  analysts  in  estimating  the  shape  parameter  of  a  system  unless  it  is  assumed 
that  shape  parameter  information  can  be  "borrowed"  from  similar  systems.  Moreover,  it  is 
seldom  worth  the  effort  needed  to  submit  the  system  to  special  testing  just  to  determine 
it's  shape  parameter.  Finally,  there  is  some  doubt  that  the  necessary  relationship  between 
expected  time  to  failure  and  the  numbers  of  failures  observed  in  a  fixed  time  exists 
between  the  exponential  and  Weibull  distributions.  Therefore,  one  area  for  further 
research  would  be  in  devising  methods  by  which  analysts  could  estimate  the  shape 
parameter  of  a  system  prior  to  testing;  see  the  "borrowing  information"  idea  suggested 
above.  Another  area  of  further  research,  which  would  be  useful  in  the  area  previously 
mentioned,  would  be  in  creating  a  database  of  system  types  and  the  associated  Weibull 
shape  parameters  of  their  inter-failure  distribution.  A  third  area  of  investigation  would  be 
in  identifying  those  combinations  of  exponential  and  Weibull  parameters  in  which  a 
suitable  relationship  exists  between  the  expected  time  to  failure  and  the  numbers  of 
failures  observed  in  a  fixed  time. 

Concerning  the  combining  of  DT  and  OT  failure  data,  the  ability  to  make  inferences 
about  the  OT  failure  rate  of  a  system  prior  to  testing  would  certainly  benefit  the  OTAs  by 
enabling  them  to  defer  testing  until  they  are  reasonably  sure  the  system  will  pass  the  test 
Moreover,  the  ability  to  augment  future  OT  failure  data  with  pseudo-data  obtained  from 
prior  DT  failure  data  would  enable  the  OTAs  to  reduce  the  amount  of  testing  resources 
needed  to  verify  reliability  thresholds.  To  date,  the  notion  of  inferring  the  relationship 
between  DT  failure  data  and  OT  failure  data  based  on  previous  similar  systems  is  untested. 
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Therefore,  further  areas  of  research  in  this  field  include  the  investigation  of  actual  data  on 
several  similar  systems  to  determine  the  feasibility  of  detecting  any  relationship  that  might 
exist  between  them  and  also  to  investigate  the  nature  of  the  relationship  which  is  found. 

Finally,  concerning  the  estimation  of  the  K  factor  using  a  hierarchical  Bayes 
methodology,  further  research  is  needed  in  order  to  develop  a  more  stable  set  of  equations 
for  the  parameters  of  the  gamma  prior.  Moreover,  one  could  consider  whether  or  not  a 
gamma  prior  is  appropriate  for  this  model,  or  whether  some  other  prior  distribution  might 
be  more  appropriate. 
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APPENDIX  A.  ALTERNATIVE  METHOD  OF  DERIVING  TIME  TO  TEST 


There  is  an  alternative  way  to  determine,  for  a  given  level  of  confidence,  the  length  of 
time,  t,  needed  to  test  a  continuous-type  system  in  order  to  demonstrate  that  the  true 
MTTF  of  the  system  is  at  least  as  large  as  the  MAV.10 

Let  the  MAV  be  denoted  by  0,  the  100a  percent  lower  confidence  limit  for  y,  the  true 
MTTF.  The  usual  test  plan  is  to  test  for  a  length  of  time,  t,  and  observe  F,  the  random 
number  of  failures  to  occur  in  (0,  t].  The  planned  test  time,  t,  is  chosen  so  that  if  one 
failure  is  observed,  the  computed  lower  confidence  limit,  yL,  will  be  equal  to  the  MAV,  0. 
The  lower  confidence  limit,  yL,  is  the  solution  for  y  in  the  equation 

F  (jVG) 

l-o  =  P(T>t)  =  P(N(t)SF)  =  £ili—  (A.1) 

It?  i! 


where  T  is  the  time  to  the  (F+l)st  failure. 

The  following  relationship  can  be  used  to  solve  for  y  in  equation  (A.l) 

Z~r= j  *-v*<i* = HxL  *  2*)  (a.2) 

j" 0  J*  1 


where  xl»  denotes  a  chi-square  random  variable  with  2a  degrees  of  freedom.  Equation 
(A.l)  then  becomes 


l0The  author  is  deeply  indebted  to  W.  Max  Woods,  Professor  of  Operations  Research,  Naval  Postgraduate 
School,  for  his  guidance  in  the  derivation  of  this  alternative  method. 


1  a  -  P  X2(F*D  -  ^ 


(A.3) 


which  implies  that 


2t_  2 

“Xa,2(F+U 


2t 

Y  =  — - 

XaJ(Ftl) 


(A.4) 


where  x^(F+1)  is  the  ath  percentile  of  the  chi-square  distribution  with  2(F+1)  degrees  of 
freedom.  Since  the  solution  for  7  is  7L,  and  the  MAY  is  0, 


9  =  Yl=-i^— 

A,a,2(F+l) 


(A.5) 


Therefore 


,  Xa.2(FH) 


(A.6) 


If  F  =  1  and  a  =  0.8,  then  Xo»a  =  5.999  and  thus,  from  equation  (A.6),  t  =  30,  as 
shown  earlier  in  section  n.A.  Note  that  equation  (A.6)  holds  for  all  positive  values  of  0 
and  F,  making  this  a  very  simple  way  to  solve  for  t  given  any  values  of  0  and  F. 
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APPENDIX  B.  CODE  FOR  RANDOM  POISSON  ARRIVAL  GENERATOR 


In  this  appendix,  the  Poisson  failure  generator  discussed  in  Chapter  IH  is  given.  The 
program  was  used  to  simulate  the  number  of  failures  during  a  given  phase  of  testing  for  a 
given  exponential  failure  rate  and  a  given  test  duration.  The  program  was  written  in 
Pascal  and  executed  using  Turbo  Pascal  6.0  (Borland,  1990),  a  commercial  Pascal 
package. 

Program  Generate_Failures; 
uses  TOOLBOX,  CRT; 

{ Author:  Timothy  P.  Anderson,  Naval  Postgraduate  School 
Date:  27  March  1994 

Purpose:  To  generate  numbers  of  failures  ihat  occur  during  a  pre  specified 
time,  TIME,  when  the  underlying  failure  distribution  is 
exponential  with  rate  LAMBDA.  This  program  is  written  in  partial 
support  of  a  thesis  by  the  author. 


var 

Lambda,  Time,  X,  CumTime,  Rand  :  real; 

I,  N,  Count :  integer; 

View:  boolean; 

OutFile :  string; 

Ch :  char; 

begin 

CRT.Clrscr; 

writeln  ('Welcome  to  Program  Generate_Failures...  Follow  instructions...'); 
writeln;  writeln; 

write  ('Specify  the  random  seed:>  ');  readln  (randseed);  writeln; 
write  ('Specify  the  exponential  failure  rate,  Lambda:>  ’); 

Lambda:=  TOOLBOX.Get_Real;  writeln; 
write  ('Specify  the  testing  time  allowed:>  ’); 

Time:=  TOOLBOX.Get_Real;  writeln; 
write  ('Specify  number  of  iterations  to  run:>  ’); 

Count:=  TOOLBOX.GetJnteger;  writeln; 

writeln  (’You  may  view  output  on  screen  or  send  it  to  a  text  file.'); 
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write  ('Do  you  want  to  view  it  on  the  screen?  (Y/N):>  '); 

View  :=  TOOLBOX.Get_Answer;  writeln;  writein;  writeln; 
if  View  =  False  then  begin 
write  ('Specify  the  output  file  (including  path):>  ’); 
read  (OutFile);  readln; 

SYSTEM.Assign  (Output,  OutFile); 
rewrite  (Output); 
end;  {if} 

for  I  :=  1  to  Count  do  begin 
X  :=  0; 

N  :=  0; 

CumTime  :=  0; 
repeat 

Rand  :=  Random; 

X  :=  -(1/Lambda)*ln(l  -  Rand); 
if  (X  +  CumTime)  >  Time  then  begin 
writeln  (N); 

CumTime  :=  CumTime  +  X; 
end 

else  begin 
N  :=  N  +  1; 

CumTime  :=  CumTime  +  X; 
end;  {if} 

until  CumTime  >  Time; 
end;  {for} 

SYSTEM.Close  (Output); 

SYSTEM. Assign  (Output,  'con'); 
rewrite  (Output); 

if  View  =  False  then  begin 

writeln  ('Output  sent  to  ',  OutFile,'...  Press  any  key  to  exit.'); 
Ch  :=  Readkey; 
end 

else  begin 

writeln  (’Press  any  key  to  exit.'); 

Ch  :=  Readkey; 
end; 

end. 
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